Phase transitions in autonomous reaction— diffusion 
systems on a one— dimensional lattice with boundaries 



Amir Aghamohammadi a Q and Mohammad Khorrami b ^ 



Department of Physics, Alzahra University, Tehran 19384, IRAN 
b Institute for Advanced Studies in Basic Sciences, 
P. 0. Box 159, Zanjan 45195, IRAN 



PACS numbers: 82.20.Mj, 02.50.Ga, 05.40.+j 
Keywords: reaction-diffusion, phase transition 

Abstract 

The family of autonomous reaction-diffusion models on a one-dimensional 
lattice with boundaries is studied. By autonomous, it is meant that the 
evolution equation for n-point functions contain only n- or less- point 
functions. It is shown that these models exhibit a static and a dynamic 
phase transition. 
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1 Introduction 



In recent years, reaction-diffusion systems have been studied by many peo- 
ple. As mean-field techniques, generally, do not give correct results for low- 
dimensional systems, people are motivated to study exactly-solvable stochastic 
models in low dimensions. Moreover, solving one-dimensional systems should in 
principle be easier. Exact results for some models on a one-dimensional lattice 
have been obtained, for example in [1-12]. Different methods have been used 
to study these models, including analytical and asymptotic methods, mean field 
methods, and large-scale numerical methods. 

Some interesting problems in non-equilibrium systems are non-equilibrium 
phase transitions described by phenomenological rate equations, and the way the 
system relaxes to its steady state. Kinetic generalizations of the Ising model, for 
example the Glauber model or the Kawasaki model, are such phenomenological 
models and have been studied extensively [13-18]. Combination of the Glauber 
and the Kawasaki dynamics has been also considered [19-22]. 

In p3l , an asymmetric generalization of the zero-temperature Glauber model 
on a lattice with boundaries was introduced. It was shown there that in the 
thermodynamic limit, the system shows two kinds of phase transitions. One of 
these is a static phase transition, the other a dynamic one. The static phase 
transition is controlled by the reaction rates, and is a discontinuous change of 
the behavior of the derivative of the stationary particle density at the end points, 
with respect to the reaction rates. The dynamic phase transition is controlled 
by the injection- and extraction- rates of the particles at the end points, and 
is a discontinuous change of the relaxation time towards the stationary config- 
uration. Other generalizations of the Glauber model consist of, for example, 
alternating-isotopic chains and alternating-bound chains (see ^4j, for exam- 
ple). People have also considered phase transitions induced through boundary 
conditions (see |25|, |2(|, for example). 

In p7f , a ten-parameter family of reaction-diffusion reactions was introduced 
for which the evolution equation of n-point functions contain only n- or less- 
point functions. What we do in this paper, is to investigate these systems on a 
finite lattice with boundaries. It will be shown that the stationary behavior of 
the system is effectively controlled by two parameters. On the one-dimensional 
boundary of this two-dimensional parameter space, there exists a phase tran- 
sition (in the thermodynamic limit, when the lattice becomes infinite), we call 
which a static phase transition. 

The relaxation time toward the stationary state of the system may depend 
on the injection- and extraction- rates at each of the boundaries. It will be 
shown that in the thermodynamic limit there are three regions, in one of them 
this time is independent of the injection- and extraction- rates, in the second 
it depends on the injection- and extraction rates at on end, and in the third it 
depends on the injection- and extraction rates at the other end. This is called 
the dynamic phase transition. 

The scheme of the paper is as follows. In section 2, autonomous reaction- 
diffusion systems with boundaries are introduced. In section 3, the static phase 
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transition of these systems is investigated. Finally, in section 4, the dynamic 
phase transition is studied. 

2 Autonomous reaction— diffusion systems with 
boundaries 

Consider a collection of particles drifting and reacting on a one-dimensional 
lattice with L sites. Each site may be occupied, A, the state corresponding to 
which is denoted by |1), or empty, 0, the state of which is denoted by |0). The 
rate of change of the state \ctj3) to the state 7<5) is H^r- It is shown in [^7| that 
the evolution equations for n-point functions are closed (involve only n- or less 
than n-point functions) iff the following conditions are satisfied by H: 

ttOI rrOO , rrOl , rrOO ttIX rrlO , rrll , rrlO n 

rrlO ttOO rrll ttOI , ttIO , r/00 , r/ll , rrOl n \ 
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one can write (P as 

u + s = u + s 

v + w = v + w. (3) 

At the end sites 1 and L, there are also injection- and extraction-rates. The 
injection- and extraction-rates in the first site are denoted by a and a', respec- 
tively. The corresponding rates in the last site are denoted by b and b 1 . It is 
then seen that 

(hk) = -(v + w + u + s)(n k ) + (v - v){n k+1 ) + (u - u)(n fc _i) 

+ w + s, 1 < k < L 
(hi) = - (v + w)(n-i) + (v — v)(n 2 ) + w + a(l - (m)) - a' (ni) 
(n L ) =-{u + s)(n L ) + (u- u)(n L -x) +s + 6(1 - (n L )) - b'(n L ), (4) 

where (n^) is the probability that the fc-th site is occupied. Comparing this with 
[ p3| , it is seen that the model considered there is a special case of this model 
with u = v — w = s = 0. 
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3 The static phase transition of the system 

The steady-state solution to (Q) is 

(n k ) = C + D 1 z* + D 2 zt L - 1 , (5) 

where z\^ satisfy 

— (it + V + W + S) + (V — v)zi t 2 + (u — U)Zi2 — Oj (6) 

and Z2 is that root the absolute value of which is greater. Let us consider this 
equation more carefully. Defining three new parameters p, q, and r through 

p :=v — v 
q :=u — u 

r :=u + s + v + w = u + s + v + w 

=u + s + v + w~u + s + v + w, (7) 

(where (||) has been used) (||) is rewritten as 

pz 2 — rz + q = 0. (8) 

Using (0) and the fact that the rates are nonnegative, it is seen that 

r > \p\,\q\,\p + i\Ap-q\- ( 9 ) 

The boundaries of the physical parameter space are thus 

r=\p + q\, and r = |p — q\. (10) 

For r — p + q, it is seen that u — v = s = w = 0, which means that AA and 00 
don't change. So, there are two equilibrium states on an infinite lattice without 
injection or extraction; cither all of the sites are occupied, or all of them are 
unoccupied. For p + q = —r, one has u — v = s — w — 0, which means that A% 
and %A don't change. So, there are two equilibrium states on an infinite lattice 
without injection or extraction; • • • • • • and • • • 0A0A 

For r = q — p, one has u = v = s = w = 0. The only nonzero rates 
are then Hq® = H±q and H®\ — Hqq. As all of the configurations can be 
converted to each other through the reactions, the equilibrium state of the 
infinite lattice without injection or extraction is unique. It is not difficult to see 
that for the special case that these four nonzero rates are equal, this state is 
• • ■ (1/2)(|0) + |1» ® (1/2)(|0) + |1)) <g> • ■ ■ . For r = p - q, one has u = v = s = 
w = 0. The only nonzero terms are then Hf® = Hqq and H\® = H®q. It is the 
same as the case r = q — p, with left and right interchanged. 

Also, if r = 0, then u = u = v = v = w = s = 0, so that (n^) = 0, 1 < 
k < L. Neglecting this trivial case, it is seen that r is positive and there are two 
parameters determining the behavior of the roots of (p|): 

P(z) := p z 2 - z + q = 0, (11) 
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where 



p := p/r 

q' := q/r. (12) 

Noting that P(l) < and P(— 1) > for r > \p + q\, it is seen that both of the 
roots of (|ll|) are real, one of them is between 1 and — 1, the other is out of this 
interval. So, \z\\ < 1 < \z%\ for r > \p + q\. In the thermodynamic limit L —> oo, 

{n k ) w C + Diz^, fc<£ 

(n fe ) pa C + D 2 Z2~ L_1 , L - fc < L. (13) 

Zi, Z2, C, Z?i, and D 2 are continuous functions of the rates. So the behavior of 
(n/c) near the ends of the lattice varies continuously with rates, and there is no 
phase transition. 

If r = p + q, one of the roots of ( pd|) is one, the other is q'/p' — q/p. If 
q > p, then Z\ = 1 and (rife) is flat for fc <C This is independent of p and g. 
However, if q < p, then Z2 — 1 and the slope of (rife) depends on the rates. For 
L — k <C i, a reverse behavior occurs. If g < p, then z 2 = 1 and (rife) is flat 
for L — k <C i. If, <? > p, then z 2 = g/p and (rife) varies with k for L — k <C i. 
To summarize, one defines two effective roots z\ and z r for sites near k — 1 (the 
left end) and k = L (the right end), respectively. We then have 

(14) 



and 




9/P, 1>P (15) 
1, q < p 



So there is a phase transition at p = q = r/2. This corresponds to u = v = s = 
w. 

If r = — p — q, one of the roots of ( pd| ) is —1 and the other is —q/p. The 
same behavior is repeated, that is 

J -1 ' -q>-p 

zi = < (16) 

{-q/p, -q<-p 

and 

z t =h 9/p ' - q> - p (17) 

[-1, -q<-p 

Again, there is a phase transition at p = q = —r/2. This corresponds to 
U — V — s — w. 
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It was seen that on two segments of the boundary of the physical parameter 
space, there exists a static phase transition. These segments (r = \p + q\) 
correspond to cases where the equilibrium state on an infinite lattice without 
injection and extraction is not unique. On the other segments of the boundary 
of the physical parameter space (r = \p — q\), where the equilibrium state is 
unique on an infinite lattice without injection and extraction, there is no static 
phase transition for the lattice with boundaries. 



4 The dynamic phase transition of the system 

The homogeneous part of (Q) can be written as 

(h k ) = h l k { ni ). (18) 

The eigenvalues and eigenvectors of the operator h satisfy 

E x k — —{v + w + u + s)x k + (v — v)x k+ i + (u - u)x k -i, k ^ 1, L 
E X\ = —{v + w + a + a')x\ + (v — v)x2, 

E x L = -(u + s + b + b')x L + {u- u)x L -u (19) 

where the eigenvalue and the eigenvector have been denoted by E and x, re- 
spectively. The solution to these is 

x k = az k 1 + I3z k 2l (20) 

where Zj's satisfy 

E = -(v + w + u + s) + (v - v)z + (u - S)z _1 , (21) 

and 

(v — v){az\ + 0z 2 ) — (E + a + a + v + w){az\ + ftz-z) = 
(u - u)^' 1 + Pz^- 1 ) -{E + b + b' + u + s)(az£ + $z%) = 0. (22) 

Defining 

5a := a + a' — (u + s), 

5b := b + b' - (v + w), (23) 

and using ( pl| ) to eliminate E, one arrives at 

[(u — u) + zi5a][(v — v)z2 +l + z^Sb] — [(u — u) + Z25a] 

x [(« - v)z[ +1 + z^5b] = 0. (24) 

This is the same as equation (15) in [p3| , with u and v replaced by u — u and 
v — v, respectively. The qualitative difference between (15) in pa] and (H) is 
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that u and v are nonnegative, whereas u — u and v — v may be negative. Defining 



v — v 



A 
B 



sgn(u — u) 



sgn(u — v) 



6 a 



4J) is simplified to 



Z 2 i+1 (l + A Zi){l + B/Z 2 ) - + A Z 2 )(l + B/Zx) = 0. (26) 



y/\(v - V)(U - u)\ 

Sb 

y/\{v - v){u - u)\' 



7L+1 



(25) 



Using (IT]); it is seen that 



z\z 2 



or 



Z X Z 2 = sgn[(u - u){v - v)}. 
The eigenvalue E is also written as 

E = —{v + w + u + s) + \J I (u — u) (v — v) [Zsgn{v — v) + Z~ 1 sgn(u — u)] 



(27) 
(28) 



= — (v + w + u + s) + sgn(w - v)^J\{u - u)(v - v)\{Zx + Z 2 ) 
= -(v + W + u + s) + sgn(u - u)y/\{u - u)(v - v)\{Z^ x + Z 2 1 ). 

Let's have a closer look at (p6|). Using (p8|) 



(29) 

is converted to a polynomial 
equation for Zj, having 2L + 2 roots. For (u — u)(v — v) > 0, Zj = 1 and 
Zj = —1 obviously satisfy (26). For (u — u)(v — v) < 0, Zj = i and Zj = —i are 



the trivial solutions of Pq). But these solutions lead to 

x k = z k {a + f3k), (30) 



not something like (20). And this form for Xk, generally does not satisfy the 
boundary conditions at fc = 1,L. So the other (nontrivial) 2L roots of ( p6| ) 
correspond to the eigenvalues of h. 

First consider the case (u — u)(v — v) > 0. If all of the roots of Zj are phases, 
then 



E < — (v + w + 11 + s) + 2\J (u — u)(v — v). 



(31) 



Equality holds if Zj = sgn(u — v) = sgn(u — u) . Normally, this is not a nontrivial 
root of d26|). But in the thermodynamic limit L — > 00, the nontrivial roots of 
26|) fill the whole unit circle. So the relaxation time for this case is 



+ V + W + S — 2 y/ (u — u)(v — v)] 1 . 



(32) 
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It is seen that it does not depend on the injection- and extraction-rates. If, 
however, some of the solutions of (E6h are not phases, then the situation is 
different. Let Z\ = Z be a root of (pq) with \Z\ > 1. In the thermodynamic 
limit L — > oo, (E26f) is turned to 



= 0, (33) 



which has the solutions 



-A,-B. (34) 



But note that we were seeking solutions with moduli greater than one. This 
shows that there is such a solution provided \A\ > 1 or \B\ > 1. If both hold, 
there exists two solutions with moduli greater than one. Suppose |^4| > 1. 
Putting Z\ = — A in (pl|), one arrives at 



E= —(v + w + u + s)— sgn(it - u)\J (u - u)(v - v)(A + A 1 ). (35) 

If sgn(u — u)A < 0, this value of E violates (|3l|), and the relaxation time is no 
more obtained from (62). In this case, 



T = [v + w + u + s + sgn(u - u)y/(u -u)(v- v)(A + A 1 )} 1 , (36) 

which is greater than (]3^), and does depend on the injection- and extraction- 
rates. This is the dynamic phase transition. The point at which this occurs is 



5a = — y{u — u)(v — v). (37) 
In terms of the injection- and extraction- rates, the transition point is 



a + a= u + s— — uj(v — v). (38) 
A similar behavior is seen at the transition point 



b + b' = v + w— \/ (u — u)(v — v). (39) 
If the injection- and extraction rates are less than this, then we have 



t = [v + w + u + s + sgn(tt - u)y/{u - u)(v - v)(B + B : )] \ (40) 

These sound quite similar to the results of p3| . But there is a difference. 
In the models studied in |23[] , either A could be less than —1 or B, and it 
was impossible that both be less than one. The reason is that there s = w = 
u — v — 0, and this means that in the physical region (where all of the rates 
are nonnegative), either the left-hand side of ( J3q) is always greater than the 
right-hand side of (|38|) , or the left-hand side of (|39[) is always greater than the 
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right-hand side of (|39|), since one of the right-hand sides is nonpositive. But it 
is not the case in the present model. Defining 



Ai := u + s — \/ {u — u)(v — v) = it + s — \/ {u — u)(v — v), 



B\ := v + w — yju — u)(v — v) — v + w — y/(u — u)(v — v), (41) 

it is seen that at least one of At or B\ are positive (apart from the special case 
it — v — s — w and u — v — s — w — 0, where both of them are zero). But it 
is also possible that both of them are positive. In general, there may be three 
phases: 

[v + w + a + a' + (u — u)(v — v)(a + a' — u — s) -1 ] -1 , region I 
[u + s + b + V + (it — u)(v — v)(b + b 1 — v — w)^ 1 }^ 1 , region II 
[v + w + u + s — 2^J (u — u)(v — S)] , otherwise 

(42) 

Regions I and II are defined as 

a + a'<Ax, a + a' - b -b' < Ax - Bi, for region I, (43) 

b + b'<Bi, a + a! -b-b' > Ax-Bi, for region II. (44) 

So the system may have two phases or three phases, depending on whether only 
one of Ai and B\ are positive or both of them are positive. In the special case 
mentioned above, the system has only one phase. This is the same Glauber 
model at zero temperature with diffusion, studied in [23[ ]. 

Next consider the case (u — u)(v — v) < 0. If all of the solutions to ( p6| ) are 
phases, then ( p9| ) shows that 

$l(E) = -(v + w + u+s), (45) 

and from that, 

t=(v + w + u+ (46) 

So in this case the relaxation time does not depend on the injection- and 
extraction- rates. Moreover, the eigenvalues of h are complex not real. This 
means that the relaxation of the system toward its stationary state is oscilla- 
tory. 

Now suppose that there exists solutions for (|2^) that are not phases. Let 
\Z\\ > 1 > \Z%\. At the thermodynamic limit, and using Z2 = —Z^ 1 , ( p6| ) is 
turned to 

(i-#)fi + #) = »• m 

The solution to this is 



ZiJ V Zi 



Zi = A, -B. (48) 
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It is obvious that to have non-phase roots, either \A\ or \B\ should be greater 
than 1. Suppose \A\ > 1. Corresponding to Z\ — A, one obtains 

E = -(v + w + u + s)- ^[(5a) 2 -\(u- u)(v - v)\]. (49) 

As \A\ > 1, the expression in the bracket is positive. If 8a < 0, then this value 
of E is greater than the right-hand side of (^) . So this value of E determines 
the relaxation time. That is, 

r ={v + w + u + s+ ((5a) _1 [((5a) 2 - \{u - u)(v - 
= [v + w + a + a' — \(u — u)(v — v)\(a + a — u — s) -1 ] -1 . (50) 



A similar argument shows that for 5b < — \/|(it — u)(v — v)\, the relaxation time 
is 

r=[u + s + b + b'- \(u-u)(v -v)\(b + b' -w -v)' 1 ]' 1 . (51) 

Finally, if both \A\ and \B\ are greater than 1, then the larger of (|50| ) and ( [5l"l) 
is the relaxation time. Defining Ai and B2 similar to (|4l|): 



Ai := u + s — \/ (u — u)(v — v) = u + s — y/(u — u)(v — v), 
B2 ■= v + w — y/(u — u)(v — v) — v + w — \J{u — u){v — v), (52) 
one arrives at 

{[v + w + a + a' + (u — u)(v — v)(a + a' — u — s) -1 ] -1 , region I 
[u + s + b + b' + (u — u)(v — v)(b + b' — v — w)^ 1 ]^ 1 , region II 
(v + w + u + s) _1 , otherwise 

(53) 



where the definitions of the regions are the same as (43) and (|44|), with Ai and 
B\ replaced by A2 and B2, respectively. Note that at least one of A2 and B2 is 
positive (apart from the special case u — v — s — w and u = v = s = w — 0, 
where both of them are zero), but it may be that both are positive. If only one 
of them is positive, the system has two phases. If both are positive, it has three 
phases. In the special case mentioned above, the system has only one phase. 

One can combine (|4^) and (^3|) in a single relation. First, note that ( pd] ) and 
(p2[) can be combined as 



A := u + s — y/\(u — u)(v — v)\ — u + s — y/\(u — u)(v — v)\, 



B := v + w — y/\(u — u)(v — v)\ = v + w — y/\(u — u)(v — v)\. (54) 

Then we can write 

[v + w + a + a' + (u — u)(v — v)(a + a' — u — s) -1 ] -1 , region I 
[u + s + b + b' + (u — u)(v — v)(b + b' — v — w)^ 1 ]^ 1 , region II 



{v + w + u + s — 2SR[y (u — u)(v — v)]} x , otherwise 

(55) 
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where the definitions of the regions are the same as ( 43 ) and (|44|) , with A\ and 
Bi replaced with A and B, respectively. 
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